#########
## Nature Human Behaviour 
## Leading Countries in Global Science Increasingly Receive More Citations than Other Countries Despite Doing Similar Research.
## https://doi.org/10.1038/s41562-022-01351-5
## Harvard Dataverse (Code and Metadata): https://doi.org/10.7910/DVN/WCOINR 
## Figures 3 through 8
## Supplemental Information Figures for Figures 3 through 8
## Data: Data_20210905
#########

#.rs.restartR()
gc(reset=T)
rm(list = ls(all.names = TRUE))

#### Libraries -----------
library("dplyr")
library("ggplot2")
library("scales")

#### Set up Parameters -----------
minimum_year <- 1980
maximum_year <- 2017
Percentile_Cutoff_Criterion_ <- 0
Citation_Window_Criterion_ <- 5
Inflation_Criterion_ <- "Field_Deflated"
Journal_Censoring_Criterion_ <- "Since_1980"
Sample_Core_Countries_ <- c("United States","China","Germany","United Kingdom","Japan","Netherlands","Switzerland")
Censored_Core_Country_List <- c("United Kingdom","Germany","Canada","Australia","New Zealand","Italy","Spain","Finland","Norway","Sweden","Austria","Portugal","France","Netherlands","Denmark","Japan","South Korea","Switzerland","Belgium","Ireland","China","Singapore","Israel","Iceland","Taiwan","Greenland","Luxembourg","Monaco","Liechtenstein","Greece")
Core_Country_List <- c("United States","United Kingdom","Germany","Canada","Australia","New Zealand","Italy","Spain","Finland","Norway","Sweden","Austria","Portugal","France","Netherlands","Denmark","Japan","South Korea","Switzerland","Belgium","Ireland","China","Singapore","Israel","Iceland","Taiwan","Greenland","Luxembourg","Monaco","Liechtenstein","Greece")

Countries_to_Drop <- c("Faroe Islands","Vatican","Jersey")

# Functions  ---------------------------------

Grand_SE <- function(x){return(sqrt(sum(x^2,na.rm=T)/length(x)^2))} #https://stats.stackexchange.com/questions/21104/calculate-average-of-a-set-numbers-with-reported-standard-errors

mean_center <- function(x){return(( (x-mean(x,na.rm=T))/sd(x,na.rm=T)))}

standard_error <- function(x){return(sd(x,na.rm=T)/sqrt(length(x)))}

# Data  ---------------------------------

MAG_Field_ID_df <- read.csv("INPUT_MAG_FieldIDs_Domain.csv") %>%
  dplyr::mutate(Domain = case_when(Domain=="Medicine" ~ "Biomedical, Behavioral,\nand Ecological Sciences",
                                   Domain=="Engineering and Computer Science" ~ "Engineering and\nComputational Sciences",
                                   Domain=="Social Sciences" ~ "Social\nSciences",
                                   Domain=="Life, Behavioral, and Ecological Sciences" ~ "Biomedical, Behavioral,\nand Ecological Sciences",
                                   Domain=="Physical and Mathematical Sciences" ~ "Physical and\nMathematical Sciences",
                                   TRUE ~ as.character(Domain))) %>%
  as.data.frame()



Regional_df <- read.csv("INPUT_R_Nation_to_Region_Core.csv") %>%
  dplyr::mutate(Region = case_when(Nation=="PUERTO_RICO" ~ "LATIN_AND_SOUTH_AMERICA",
                                   Nation=="CUBA" ~ "LATIN_AND_SOUTH_AMERICA",
                                   Nation=="GREENLAND" ~ "WESTERN_EUROPE",
                                   Nation=="FRENCH_POLYNESIA" ~ "OCEANIA",
                                   Nation=="BRUNEI_DARUSSALAM" ~ "SOUTH_EAST_ASIA",
                                   Nation=="COMOROS" ~ "SUB_AFRICA",
                                   Nation=="CABO_VERDE" ~ "SUB_AFRICA",
                                   Nation=="KIRIBATI" ~ "OCEANIA",
                                   Nation=="LAO_PDR" ~ "SOUTH_EAST_ASIA",
                                   Nation=="MOLDOVA" ~ "EASTERN_EUROPE_USSR",
                                   Nation=="MAURITANIA" ~ "MENA",
                                   Nation=="NEW_CALEDONIA" ~ "OCEANIA",
                                   Nation=="SAN_MARINO" ~ "WESTERN_EUROPE",
                                   Nation=="SAO_TOME_AND_PRINCIPE" ~ "SUB_AFRICA",
                                   Nation=="TUVALU" ~ "OCEANIA",
                                   Nation=="ST_VINCENT_AND_THE_GRENADINES" ~ "CARIBBEAN",
                                   Nation=="MAURITIUS" ~ "SUB_AFRICA",
                                   Nation=="RUSSIA" ~ "CENTRAL_ASIA",
                                   Nation=="KAZAKHSTAN" ~ "CENTRAL_ASIA",
                                   TRUE ~ as.character(Region))) %>%
  dplyr::mutate(Region = case_when(
    Region=="CARIBBEAN" ~ "Latin America\nand Caribbean",
    Region=="WESTERN_EUROPE" ~ "Europe",
    Region=="SUB_AFRICA" ~ "MENA and\nSub-Saharan Africa",
    Region=="EASTERN_EUROPE_USSR" ~ "Europe",
    Region=="EAST_ASIA" ~ "Asia",
    Region=="CENTRAL_ASIA" ~ "Asia",
    Region=="SOUTH_EAST_ASIA" ~ "Asia",
    Region=="LATIN_AND_SOUTH_AMERICA" ~ "Latin America\nand Caribbean",
    Region=="NORTH_AMERICA" ~ "US and Canada",
    Region=="MENA" ~ "MENA and\nSub-Saharan Africa",
    Region=="OCEANIA" ~ "Oceania"
  ))%>%
  dplyr::mutate(Country = stringr::str_to_title(stringr::str_replace(as.character(Nation),"_"," "))) %>%
  unique() %>%
  dplyr::mutate(Country = case_when(Country=="Republic Of_the_congo"~"Republic of the Congo",
                                    Country=="Saint Vincent_and_the_grenadines"~"Saint Vincent",
                                    Country=="State Of_palestine"~"Palestine",
                                    Country=="Democratic Republic_of_the_congo"~"Democractic Republic of the Congo",
                                    Country=="Papua New_guinea"~"Papua New Guinea",
                                    Country=="Trinidad And_tobago"~"Trinidad and Tobago",
                                    Country=="Saint Kitts_and_nevis"~"Saint Kitts and Nevis",
                                    Country=="United Arab_emirates"~"United Arab Emirates",
                                    Country=="Antigua And_barbuda"~"Antigua and Barbuda",
                                    Country=="British Virgin_islands"~"British Virgin Islands",
                                    Country=="Bosnia And_herzegovina"~"Bosnia and Herzegovina",
                                 TRUE~as.character(Country)))%>%
  as.data.frame() %>%
  select("Country","Region")

# Plot 3 | Delta ---------------------------------

Delta_Degree_Centrality_df <- plyr::ldply(list.files(pattern="OUTPUT_Python_Field_Delta_Degree_Centrality_KLD_Citation_Corpus_RAKE_and_GoogleAPI_EnglishOnly_*",
                                                     full.names=T), 
                                          read.csv) %>%  
  dplyr::mutate(Country = stringr::str_to_title(stringr::str_replace(as.character(Country),"_"," "))) %>%
  dplyr::mutate(Country = case_when(Country=="Republic Of_the_congo"~"Republic of the Congo",
                                    Country=="Saint Vincent_and_the_grenadines"~"Saint Vincent",
                                    Country=="State Of_palestine"~"Palestine",
                                    Country=="Democratic Republic_of_the_congo"~"Democractic Republic of the Congo",
                                    Country=="Papua New_guinea"~"Papua New Guinea",
                                    Country=="Trinidad And_tobago"~"Trinidad and Tobago",
                                    Country=="Saint Kitts_and_nevis"~"Saint Kitts and Nevis",
                                    Country=="United Arab_emirates"~"United Arab Emirates",
                                    Country=="Antigua And_barbuda"~"Antigua and Barbuda",
                                    Country=="British Virgin_islands"~"British Virgin Islands",
                                    Country=="Bosnia And_herzegovina"~"Bosnia and Herzegovina",
                                    Country=="Czechia"~"Czech Republic",
                                    TRUE~as.character(Country)))%>%
  as.data.frame()%>%
  subset(!(Country%in%Countries_to_Drop))

# Table Names of Fields ---------------------------------

Field_Tables_df <- Delta_Degree_Centrality_df %>%
  dplyr::select(Discipline) %>%
  unique() %>%
  merge(.,MAG_Field_ID_df %>% dplyr::select(FieldID,Domain,Name) %>% dplyr::rename(Discipline=FieldID),
        by=c("Discipline"),all.x=T) %>%
  dplyr::select(Domain,Name) %>%
  unique() %>%
  dplyr::mutate(Name = case_when(Name=="Human‚Äìcomputer interaction"~"Human-computer interaction",
                                 TRUE~as.character(Name)))%>%
  dplyr::mutate(Name = tools::toTitleCase(as.character(Name))) %>%
  dplyr::group_by(Domain) %>%
  arrange(Name)%>%
  dplyr::summarise(Name = paste0(Name, collapse = " | ")) %>%
  as.data.frame()

#write.csv(Field_Tables_df,"TABLE_R_Field_Names_by_Domain.csv",row.names = F)

# Table List of Countries (Region) ---------------------------------

# Nation_Tables_df_Region <- KLD_Degree_Centrality_df %>%
#   dplyr::select(Country) %>%
#   unique()%>%
#   merge(.,Regional_df%>%
#           dplyr::rename(Country=Nation)%>%
#           dplyr::mutate(Country = tools::toTitleCase(tolower(as.character(Country)))),
#         by=c("Country"),all.x=T) %>%
#   unique() %>%
#   dplyr::mutate(Name = tools::toTitleCase(tolower(as.character(Country)))) %>%
#   dplyr::group_by(Region) %>%
#   dplyr::summarise(Name = paste0(Name, collapse = " | ")) %>%
#   as.data.frame()
# 
# write.csv(Nation_Tables_df_Region,"TABLE_R_Country_Names_by_Region.csv",row.names = F)

# Table List of Countries (Core-Periphery) ---------------------------------

Nation_Tables_df_CorePeriphery <- Delta_Degree_Centrality_df %>%
  dplyr::select(Country) %>%
  unique()%>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Full_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Core","Periphery")))%>%
  dplyr::mutate(Name = tools::toTitleCase(tolower(as.character(Country)))) %>%
  dplyr::mutate(Name = case_when(Name=="Republic Of_the_congo"~"Republic of the Congo",
                                 Name=="Saint Vincent_and_the_grenadines"~"Saint Vincent",
                                 Name=="State Of_palestine"~"Palestine",
                                 Name=="Democratic Republic_of_the_congo"~"Democractic Republic of the Congo",
                                 Name=="Papua New_guinea"~"Papua New Guinea",
                                 Name=="Trinidad And_tobago"~"Trinidad and Tobago",
                                 Name=="Saint Kitts_and_nevis"~"Saint Kitts and Nevis",
                                 Name=="United Arab_emirates"~"United Arab Emirates",
                                 Name=="Antigua And_barbuda"~"Antigua and Barbuda",
                                 Name=="British Virgin_islands"~"British Virgin Islands",
                                 Name=="Bosnia And_herzegovina"~"Bosnia and Herzegovina",
                                 TRUE~as.character(Name)))%>%
  dplyr::group_by(Core_Periphery) %>%
  arrange(Name)%>%
  dplyr::summarise(Name = paste0(Name, collapse = " | ")) %>%
  as.data.frame()

#write.csv(Nation_Tables_df_CorePeriphery,"TABLE_R_Country_Names_by_CorePeriphery.csv",row.names = F)

# Plot 3A | Delta | Data and Plot ---------------------------------

Delta_Degree_Centrality_df_Trends <- Delta_Degree_Centrality_df %>%
  dplyr::group_by(Country,Year,Percentile_Cutoff,Citation_Type,Journal_Censoring) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE=sd(InDelta_Citation_ij,na.rm=T)/sqrt(length(InDelta_Citation_ij))) %>%
  as.data.frame()%>%
  subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
  tidyr::drop_na()%>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile"))

Delta_Degree_Centrality_df_Trends_Main <- Delta_Degree_Centrality_df_Trends %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)

Delta_Degree_Centrality_df_Trends_Main_Core_Periphery <- Delta_Degree_Centrality_df %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_) %>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Core_Country_List ~ "Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::group_by(Core_Periphery,Year)%>%
  dplyr::summarise(InDelta_Citation_ij_Mean = mean(InDelta_Citation_ij,na.rm=T),
                   InDelta_Citation_ij_SE = sd(InDelta_Citation_ij,na.rm=T)/sqrt(length(InDelta_Citation_ij)))%>%
  as.data.frame()


plot_Delta_over_CountrySpecific <- ggplot(data=subset(Delta_Degree_Centrality_df_Trends_Main,
                                                      Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_),
                                          aes(x=Year,
                                              y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  xlim(c(1980,2025))+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Trends_Main,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Trends_Main,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Trends_Main,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  fill=Country,
                  colour=Country),               
                  alpha=0.2)+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_Main,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2025),
                            #size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_Main,
  #                                       Year==2010 &  Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country),
  #                           fontface = 'bold', color = 'white',segment.color = "black",show.legend = F)+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()

# Plot 3C | Delta Domain | Data and Plot ---------------------------------

Delta_Degree_Centrality_df_Domain <- Delta_Degree_Centrality_df %>%
  merge(.,MAG_Field_ID_df %>%
          dplyr::rename(Discipline=FieldID) %>%
          select(Discipline,Domain),
        by=c("Discipline"),
        all.x=T) %>%
  dplyr::group_by(Domain,Country,Year,Percentile_Cutoff,Citation_Window,Citation_Type,Journal_Censoring) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE=sd(InDelta_Citation_ji,na.rm=T)/sqrt(length(InDelta_Citation_ji))) %>%
  as.data.frame() %>%
  subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
  tidyr::drop_na()%>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile"))

Delta_Degree_Centrality_df_Domain_Main <- Delta_Degree_Centrality_df_Domain %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_ &
           Journal_Censoring == Journal_Censoring_Criterion_ & 
           Citation_Window==Citation_Window_Criterion_ & 
           Citation_Type==Inflation_Criterion_ & 
           Year>=minimum_year & 
           Year<=maximum_year-Citation_Window_Criterion_)


plot_Delta_over_CountrySpecific_Domain <- ggplot(data=subset(Delta_Degree_Centrality_df_Domain_Main,
                                                             Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_),
                                                 aes(x=Year,
                                                     y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  
  xlim(c(1980,2035))+
  #ylim(c(-1.0,3))+ #c(-1,2.74)
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Domain_Main,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Domain_Main,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Domain_Main,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,fill=Country),               
                  alpha=0.2)+
  
   ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Domain_Main,
                                         Year==2012 &  Country%in%Sample_Core_Countries_),
                             aes(fill=Country,label=Country),
                             #hjust = -0.05,
                             xlim = c(1980,2035),
                             size = 3,
                             force= 1,
                             nudge_x = 30,
                             direction = 'y',
                             fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(~Domain,scale="free_y")


# Plot 3B | Delta Core and Periphery | Data and Plot ---------------------------------

Delta_Degree_Centrality_df_Core <- Delta_Degree_Centrality_df %>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))%>%
  dplyr::group_by(Core_Periphery,Year,Percentile_Cutoff,Citation_Window,Citation_Type,Journal_Censoring) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE=sd(InDelta_Citation_ji,na.rm=T)/sqrt(length(InDelta_Citation_ji))) %>%
  as.data.frame() %>%
  subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
  tidyr::drop_na()%>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile")) %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_ &
           Journal_Censoring == Journal_Censoring_Criterion_ & 
           Citation_Window==Citation_Window_Criterion_ & 
           Citation_Type==Inflation_Criterion_ & 
           Year>=minimum_year & 
           Year<=maximum_year-Citation_Window_Criterion_)

plot_Delta_CorePeriphery <- ggplot(data=Delta_Degree_Centrality_df_Core,
                                   aes(x=Year,y=InDelta_Citation_ij_Mean,group=Core_Periphery))+
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))  +
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  ggtitle("")+
  geom_ribbon( aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                   ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                   colour=Core_Periphery,
                   fill=Core_Periphery),
               alpha=0.2)+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Core,
  #                                        Year==2012),
  #                            aes(fill=Core_Periphery,label=Core_Periphery,size=10), 
  #                            fontface = 'bold', color = 'white',segment.color = "black")+
  theme(plot.title = element_text(hjust = 0.5))


# Plot 3D | Delta Core and Periphery Domain | Data and Plot ---------------------------------

Delta_Degree_Centrality_df_Core_Domain <- Delta_Degree_Centrality_df %>%
  merge(.,MAG_Field_ID_df %>%
          dplyr::rename(Discipline=FieldID) %>%
          select(Discipline,Domain),
        by=c("Discipline"),
        all.x=T) %>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  subset(Core_Periphery!="Drop")%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))%>%
  dplyr::group_by(Domain,Core_Periphery,Year,Percentile_Cutoff,Citation_Window,Citation_Type,Journal_Censoring) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE=sd(InDelta_Citation_ij,na.rm=T)/sqrt(length(InDelta_Citation_ij))) %>%
  as.data.frame() %>%
  subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
  tidyr::drop_na()%>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile")) %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_ &
           Journal_Censoring == Journal_Censoring_Criterion_ & 
           Citation_Window==Citation_Window_Criterion_ & 
           Citation_Type==Inflation_Criterion_ & 
           Year>=minimum_year & 
           Year<=maximum_year-Citation_Window_Criterion_)

plot_Delta_CorePeriphery_Domain <- ggplot(data=Delta_Degree_Centrality_df_Core_Domain,
                                          aes(x=Year,y=InDelta_Citation_ij_Mean,group=Core_Periphery))+
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))  +
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  ggtitle("")+
  geom_ribbon( aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                   ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                   colour=Core_Periphery,
                   fill=Core_Periphery),
               alpha=0.2)+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Core_Domain,
  #                                        Year==2012),
  #                            aes(fill=Core_Periphery,label=Core_Periphery,size=5), 
  #                            fontface = 'bold', color = 'white',segment.color = "black")+
  theme(plot.title = element_text(hjust = 0.5))+
  facet_grid(~Domain)

# Core and Periphery First Apperance Impact on Delta | Data (and Domain Data) ---------------------------

FirstApperance_Decade_df <- Delta_Degree_Centrality_df %>%
  select(Discipline,Country,Year) %>%
  unique()%>%
  dplyr::group_by(Discipline,Country)%>%
  dplyr::summarise(Year = min(Year, na.rm = T))%>%
  dplyr::mutate(FirstApperance_Decade = case_when(Year>=1980 & Year<1985 ~ "1980 to 1984",
                                              Year>=1985 & Year<1990 ~ "1985 to 1989",
                                              Year>=1990 & Year<1995 ~ "1990 to 1994",
                                              Year>=1995 & Year<2000 ~ "1995 to 1999",
                                              Year>=2000 & Year<2005 ~ "2000 to 2004",
                                              Year>=2005 & Year<=2012 ~ "2005 to 2012",
                                              TRUE ~ as.character(NA)))%>%
  as.data.frame()%>%
  tidyr::drop_na()%>%
  
  select(Discipline,Country,FirstApperance_Decade)

Delta_Degree_Centrality_df_Core_by_Decade_FirstApperance <- Delta_Degree_Centrality_df_Core_by_Decade_DOMAIN_FirstApperance <- data.frame()

for( i in c(1,2,3,4,5,6)){

  if(i==1){decade_list <- c("1980 to 1984")
  decade_title <- c("First Apperance\n1980 to 1984")}
  if(i==2){decade_list <- c("1980 to 1984","1985 to 1989")
  decade_title <- c("First Apperance\n1980 to 1989")}
  if(i==3){decade_list <- c("1980 to 1984","1985 to 1989","1990 to 1994")
  decade_title <- c("First Apperance\n1980 to 1994")}
  if(i==4){decade_list <- c("1980 to 1984","1985 to 1989","1990 to 1994","1995 to 1999")
  decade_title <- c("First Apperance\n1980 to 1999")}
  if(i==5){decade_list <- c("1980 to 1984","1985 to 1989","1990 to 1994","1995 to 1999","2000 to 2004")
  decade_title <- c("First Apperance\n1980 to 2004")}
  if(i==6){decade_list <- c("1980 to 1984","1985 to 1989","1990 to 1994","1995 to 1999","2000 to 2004","2005 to 2012")
  decade_title <- c("First Apperance\n1980 to 2012")}

  Delta_Degree_Centrality_df_Core_by_Decade_FirstApperance_Individual <- Delta_Degree_Centrality_df %>%
  merge(.,FirstApperance_Decade_df,
        by=c("Discipline","Country"),all.x=T)%>%
  subset(FirstApperance_Decade%in%c(decade_list))%>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
    
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))%>%
  dplyr::group_by(Core_Periphery,Year,Percentile_Cutoff,Citation_Window,Citation_Type,Journal_Censoring) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE=sd(InDelta_Citation_ji,na.rm=T)/sqrt(length(InDelta_Citation_ji))) %>%
  as.data.frame() %>%
  subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
  tidyr::drop_na()%>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile")) %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_ &
           Journal_Censoring == Journal_Censoring_Criterion_ & 
           Citation_Window==Citation_Window_Criterion_ & 
           Citation_Type==Inflation_Criterion_ & 
           Year>=minimum_year & 
           Year<=maximum_year-Citation_Window_Criterion_) %>%
  dplyr::mutate(Decade_Title = decade_title)
  
  Delta_Degree_Centrality_df_Core_by_Decade_DOMAIN_FirstApperance_Individual <- Delta_Degree_Centrality_df %>%
    merge(.,FirstApperance_Decade_df,
          by=c("Discipline","Country"),all.x=T)%>%
    subset(FirstApperance_Decade%in%c(decade_list))%>%
    merge(.,MAG_Field_ID_df %>%
            dplyr::rename(Discipline=FieldID) %>%
            select(Discipline,Domain),
          by=c("Discipline"),
          all.x=T) %>%
    dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                             TRUE ~ "Periphery"))%>%
    dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))%>%
    dplyr::group_by(Domain,Core_Periphery,Year,Percentile_Cutoff,Citation_Window,Citation_Type,Journal_Censoring) %>%
    dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                     InDelta_Citation_ij_SE=sd(InDelta_Citation_ji,na.rm=T)/sqrt(length(InDelta_Citation_ji))) %>%
    as.data.frame() %>%
    subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
    tidyr::drop_na()%>%
    dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                       Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                       Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                       Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile")) %>%
    subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_ &
             Journal_Censoring == Journal_Censoring_Criterion_ & 
             Citation_Window==Citation_Window_Criterion_ & 
             Citation_Type==Inflation_Criterion_ & 
             Year>=minimum_year & 
             Year<=maximum_year-Citation_Window_Criterion_) %>%
    dplyr::mutate(Decade_Title = decade_title)  

  Delta_Degree_Centrality_df_Core_by_Decade_FirstApperance <- rbind(Delta_Degree_Centrality_df_Core_by_Decade_FirstApperance,
                                                                    Delta_Degree_Centrality_df_Core_by_Decade_FirstApperance_Individual)

  Delta_Degree_Centrality_df_Core_by_Decade_DOMAIN_FirstApperance <- rbind(Delta_Degree_Centrality_df_Core_by_Decade_DOMAIN_FirstApperance,
                                                                    Delta_Degree_Centrality_df_Core_by_Decade_DOMAIN_FirstApperance_Individual)
  
}

# Core and Periphery First Apperance Impact on Delta | Plot ---------------------------

SM_plot_CorePeriphery_FirstApperance_Impact_on_Delta <- ggplot(data=Delta_Degree_Centrality_df_Core_by_Decade_FirstApperance,
                                                             aes(x=Year,y=InDelta_Citation_ij_Mean,group=Core_Periphery))+
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))  +
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  ggtitle("")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  geom_ribbon( aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                   ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                   colour=Core_Periphery,
                   fill=Core_Periphery),
               alpha=0.2)+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Core_by_Decade_FirstApperance,
  #                                       Year==2012),
  #                           aes(fill=Core_Periphery,label=Core_Periphery,size=10), 
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  theme(plot.title = element_text(hjust = 0.5))+
  facet_wrap(~Decade_Title)

# Core and Periphery First Apperance Impact on Delta | Domain Plot ---------------------------

SM_plot_CorePeriphery_FirstApperance_Impact_on_Delta_Domain <- ggplot(data=Delta_Degree_Centrality_df_Core_by_Decade_DOMAIN_FirstApperance,
                                                               aes(x=Year,y=InDelta_Citation_ij_Mean,group=Core_Periphery))+
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))  +
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  ggtitle("")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  geom_ribbon( aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                   ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                   colour=Core_Periphery,
                   fill=Core_Periphery),
               alpha=0.2)+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Core_by_Decade_DOMAIN_FirstApperance,
  #                                       Year==2012),
  #                           aes(fill=Core_Periphery,label=Core_Periphery,size=10), 
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  theme(plot.title = element_text(hjust = 0.5))+
  facet_grid(Domain~Decade_Title)


# Percentile Cutoff (Delta) -----------------
# Percentile Cutoff (Delta) | Data and Plot ---------------------------------

Delta_Degree_Centrality_df_Trends_SM <- Delta_Degree_Centrality_df_Trends %>%
  subset(Percentile_Cutoff!=Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)

SM_Cutoff_plot_Delta_over_CountrySpecific <- ggplot(data=Delta_Degree_Centrality_df_Trends_SM,
                                                    aes(x=Year,
                                                        y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Trends_SM,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Trends_SM,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Trends_SM,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,fill=Country),alpha=0.2)+
  
  xlim(c(1980,2025))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_SM,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2025),
                            #size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_SM,
  #                                       Year==2010 &  Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(~Percentile_Cutoff_Labels)



# Percentile Cutoff (Delta) | Domain Data and Plot ---------------------------------

Delta_Degree_Centrality_df_Domain_SM <- Delta_Degree_Centrality_df_Domain %>%
  subset(Percentile_Cutoff!=Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)

SM_Cutoff_plot_Delta_over_CountrySpecific_Domain <- ggplot(data=Delta_Degree_Centrality_df_Domain_SM,
                                                           aes(x=Year,
                                                               y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  #ylim(c(-1.0,2.75))+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Domain_SM,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Domain_SM,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Domain_SM,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,fill=Country),alpha=0.2)+
  
  xlim(c(1980,2035))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Domain_SM,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2035),
                            size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Domain_SM,
  #                                       Year==2010 &  Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Percentile_Cutoff_Labels~Domain,scale="free_y")

# Deflated vs. Standard (Delta) -------------------------
# Deflated vs. Standard (Delta) | Data and Plot -------------------------

Delta_Degree_Centrality_df_Trends_Deflated <- Delta_Degree_Centrality_df_Trends %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  #subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Citation_Type_Label = case_when(Citation_Type=="No_Censoring"~"No Deflation",
                                                Citation_Type=="Field_Deflated"~"Field Deflated",
                                                Citation_Type=="Dyad_Deflated"~"Receiver Deflated"))%>%
  as.data.frame()

SM_Deflated_plot_Delta_over_CountrySpecific <- ggplot(data=Delta_Degree_Centrality_df_Trends_Deflated,
                                                    aes(x=Year,
                                                        y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Trends_Deflated,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Trends_Deflated,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Trends_Deflated,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,fill=Country),alpha=0.2)+
  
  xlim(c(1980,2025))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_Deflated,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2025),
                            #size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_Deflated,
  #                                       Year==2010 &  Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(~Citation_Type_Label)

# Deflated vs. Standard (Delta) | Domain Data and Plot -------------------------

Delta_Degree_Centrality_df_Trends_Deflated_DOMAIN <- Delta_Degree_Centrality_df_Domain %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  #subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Citation_Type_Label = case_when(Citation_Type=="No_Censoring"~"No Deflation",
                                                Citation_Type=="Field_Deflated"~"Field Deflated",
                                                Citation_Type=="Dyad_Deflated"~"Receiver Deflated"))%>%
  as.data.frame()

SM_Deflated_plot_Delta_over_CountrySpecific_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_Trends_Deflated_DOMAIN,
                                                      aes(x=Year,
                                                          y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Trends_Deflated_DOMAIN,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Trends_Deflated_DOMAIN,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Trends_Deflated_DOMAIN,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,fill=Country),alpha=0.2)+
  
  xlim(c(1980,2035))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_Deflated_DOMAIN,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2035),
                            size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_Deflated_DOMAIN,
  #                                       Year==2010 &  Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Domain~Citation_Type_Label)

# English vs. Non-English (Delta) -------------------------

Delta_Degree_Centrality_df_All <- plyr::ldply(list.files(pattern="OUTPUT_Python_Field_Delta_Degree_Centrality_KLD_Citation_Corpus_RAKE_and_GoogleAPI_All_*",
                                                         full.names=T), 
                                              read.csv)  %>%
  dplyr::mutate(Country = stringr::str_to_title(stringr::str_replace(as.character(Country),"_"," "))) %>%
  subset(Percentile_Cutoff==0)%>%
  subset(Citation_Window==5)%>%
  subset(Citation_Type=="Deflated")%>%
  dplyr::select(Country,Year,Discipline,OutDelta_Citation_ji)%>%
  dplyr::rename(InDelta_Citation_ij = OutDelta_Citation_ji) %>%
  dplyr::mutate(Language = "All")%>%
  rbind(.,Delta_Degree_Centrality_df%>%
          subset(Citation_Type=="Field_Deflated")%>%
          subset(Citation_Window==5)%>%
          subset(Journal_Censoring=="Since_1980")%>%
          subset(Percentile_Cutoff==0)%>%
          dplyr::select(Country,Year,Discipline,InDelta_Citation_ij)%>%
          dplyr::mutate(Language = "English_Only"))%>%
  as.data.frame()%>%
  merge(.,subset(MAG_Field_ID_df,select=c("FieldID","Domain")),by.x=c("Discipline"),by.y=c("FieldID"),how="left")%>%
  subset(Year>=minimum_year & Year<=(maximum_year-Citation_Window_Criterion_)) %>% 
  dplyr::mutate(Language = case_when(as.character(Language)=="English_Only"~"English Only",
                                     TRUE~"English and\nTranslated"))


# English vs. Non-English (Delta) | Data and Plot -------------------------

Delta_Degree_Centrality_df_Language <- Delta_Degree_Centrality_df_All %>%
  dplyr::group_by(Year,Discipline,Language) %>%
  dplyr::mutate(InDelta_Citation_ij = mean_center(InDelta_Citation_ij)) %>%
  dplyr::group_by(Year,Country,Language) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE = sd(InDelta_Citation_ij,na.rm=T)/sqrt(length(InDelta_Citation_ij)))

SM_Language_plot_Delta_over_CountrySpecific <- ggplot(data=Delta_Degree_Centrality_df_Language ,
                                                      aes(x=Year,
                                                          y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  ### Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Language,
                         Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Language,
                        Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Language,
                          Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,
                  fill=Country),alpha=0.2)+
  
  xlim(c(1980,2025))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Language,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2025),
                            #size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Language,
  #                                       Year==2010 & Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country), fontface = 'bold', color = 'white',segment.color = "black")+
  
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_wrap(~Language)

# English vs. Non-English (Delta) | Domain Data and Plot -------------------------

Delta_Degree_Centrality_df_DOMAIN_Language <- Delta_Degree_Centrality_df_All %>%
  dplyr::group_by(Year,Discipline,Language,Domain) %>%
  dplyr::mutate(InDelta_Citation_ij = mean_center(InDelta_Citation_ij)) %>%
  dplyr::group_by(Year,Country,Language,Domain) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE = sd(InDelta_Citation_ij,na.rm=T)/sqrt(length(InDelta_Citation_ij)))

SM_Language_plot_Delta_over_CountrySpecific_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_DOMAIN_Language ,
                                                             aes(x=Year,
                                                                 y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  ### Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_DOMAIN_Language,
                         Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_DOMAIN_Language,
                        Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_DOMAIN_Language,
                          Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,
                  fill=Country),
              alpha=0.2)+
  
  xlim(c(1980,2035))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_DOMAIN_Language,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2035),
                            size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_DOMAIN_Language,
  #                                       Year==2010 & Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country), fontface = 'bold', color = 'white',segment.color = "black")+
  
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Domain~Language)

# Journal Censored vs. All (Delta) -------------------------
# Journal Censored vs. All (Delta) | Data and Plot -------------------------

Delta_Degree_Centrality_df_Trends_JournalCensored <- Delta_Degree_Centrality_df_Trends %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  #subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Journal_Censoring_Label = case_when(Journal_Censoring=="No_Censoring"~"No Journal\nCensoring",
                                                    TRUE~"Journals Censored\nSince 1980"))%>%
  as.data.frame()

SM_JournalCensored_plot_Delta_over_CountrySpecific <- ggplot(data=Delta_Degree_Centrality_df_Trends_JournalCensored,
                                                             aes(x=Year,
                                                                 y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,fill=Country),alpha=0.2)+
  
  xlim(c(1980,2025))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2025),
                            #size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored,
  #                                       Year==2010 &  Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(~Journal_Censoring_Label)

# Journal Censored vs. All (Delta) | Domain Data and Plot -------------------------

Delta_Degree_Centrality_df_Trends_JournalCensored_DOMAIN <- Delta_Degree_Centrality_df_Domain %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  #subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Journal_Censoring_Label = case_when(Journal_Censoring=="No_Censoring"~"No Journal\nCensoring",
                                                    TRUE~"Journals Censored\nSince 1980"))%>%
  as.data.frame()

SM_JournalCensored_plot_Delta_over_CountrySpecific_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_Trends_JournalCensored_DOMAIN,
                                                                    aes(x=Year,
                                                                        y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored_DOMAIN,
                         Year>=minimum_year& Year<=maximum_year-Citation_Window_Criterion_ &  Country%in%Sample_Core_Countries_),
             aes(colour=Country))+
  geom_line(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored_DOMAIN,
                        Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
            aes(colour=Country))+
  geom_ribbon(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored_DOMAIN,
                          Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_ & Country%in%Sample_Core_Countries_),
              aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Country,fill=Country),alpha=0.2)+
  
  xlim(c(1980,2035))+
  ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored_DOMAIN,
                                        Year==2012 &  Country%in%Sample_Core_Countries_),
                            aes(fill=Country,label=Country),
                            #hjust = -0.05,
                            xlim = c(1980,2035),
                            size = 3,
                            force= 1,
                            nudge_x = 30,
                            direction = 'y',
                            fontface = 'bold', color ='white',segment.color = "black",show.legend = F)+
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_Trends_JournalCensored_DOMAIN,
  #                                       Year==2010 &  Country%in%Sample_Core_Countries_),
  #                           aes(fill=Country,label=Country),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Domain~Journal_Censoring_Label)

# SM (Delta - Core and Periphery) -----------------

Delta_Degree_Centrality_df_CorePeriphery_SM <- Delta_Degree_Centrality_df %>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))%>%
  dplyr::group_by(Core_Periphery,Year,Percentile_Cutoff,Citation_Window,Citation_Type,Journal_Censoring) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE=sd(InDelta_Citation_ji,na.rm=T)/sqrt(length(InDelta_Citation_ji))) %>%
  as.data.frame() %>%
  subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
  tidyr::drop_na()%>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile"))

Delta_Degree_Centrality_df_CorePeriphery_DOMAIN_SM <- Delta_Degree_Centrality_df %>%
  merge(.,MAG_Field_ID_df %>%
          dplyr::rename(Discipline=FieldID) %>%
          select(Discipline,Domain),
        by=c("Discipline"),
        all.x=T) %>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))%>%
  dplyr::group_by(Domain,Core_Periphery,Year,Percentile_Cutoff,Citation_Window,Citation_Type,Journal_Censoring) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE=sd(InDelta_Citation_ji,na.rm=T)/sqrt(length(InDelta_Citation_ji))) %>%
  as.data.frame() %>%
  subset(Percentile_Cutoff %in% c(0,0.25,0.5,0.75)) %>%
  tidyr::drop_na()%>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.5 ~ "Above 50th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile"))

# Percentile Cutoff (Delta - Core and Periphery) -----------------
# Percentile Cutoff (Delta - Core and Periphery) | Data and Plot ---------------------------------


Delta_Degree_Centrality_df_CorePeriphery_SM_Percentile <- Delta_Degree_Centrality_df_CorePeriphery_SM %>%
  subset(Percentile_Cutoff!=Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)

SM_Cutoff_plot_Delta_over_CorePeriphery <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_SM_Percentile,
                                                  aes(x=Year,
                                                      y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,fill=Core_Periphery),alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_SM_Percentile,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(~Percentile_Cutoff_Labels)



# Percentile Cutoff (Delta - Core and Periphery) | Domain Data and Plot ---------------------------------

Delta_Degree_Centrality_df_CorePeriphery_Percentile_DOMAIN_SM <- Delta_Degree_Centrality_df_CorePeriphery_DOMAIN_SM %>%
  subset(Percentile_Cutoff!=Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)

SM_Cutoff_plot_Delta_over_CorePeriphery_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_Percentile_DOMAIN_SM,
                                                         aes(x=Year,
                                                             y=InDelta_Citation_ij_Mean))+
  
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  #ylim(c(-1.0,2.75))+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  # Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,fill=Core_Periphery),alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_Percentile_DOMAIN_SM,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Percentile_Cutoff_Labels~Domain,scale="free_y")

# Deflated vs. Standard (Delta - Core and Periphery) -------------------------
# Deflated vs. Standard (Delta - Core and Periphery) | Data and Plot -------------------------

Delta_Degree_Centrality_df_CorePeriphery_Deflated <- Delta_Degree_Centrality_df_CorePeriphery_SM %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  #subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Citation_Type_Label = case_when(Citation_Type=="No_Censoring"~"No Deflation",
                                                Citation_Type=="Field_Deflated"~"Field Deflated",
                                                Citation_Type=="Dyad_Deflated"~"Receiver Deflated"))%>%
  as.data.frame()

SM_Deflated_plot_Delta_over_CorePeriphery <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_Deflated,
                                                    aes(x=Year,
                                                        y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,fill=Core_Periphery),alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_Deflated,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(~Citation_Type_Label)

# Deflated vs. Standard (Delta - Core and Periphery) | Domain Data and Plot -------------------------

Delta_Degree_Centrality_df_CorePeriphery_Deflated_DOMAIN <- Delta_Degree_Centrality_df_CorePeriphery_DOMAIN_SM %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  #subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Citation_Type_Label = case_when(Citation_Type=="No_Censoring"~"No Deflation",
                                                Citation_Type=="Field_Deflated"~"Field Deflated",
                                                Citation_Type=="Dyad_Deflated"~"Receiver Deflated"))%>%
  as.data.frame()

SM_Deflated_plot_Delta_over_CorePeriphery_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_Deflated_DOMAIN,
                                                           aes(x=Year,
                                                               y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,fill=Core_Periphery),alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_Deflated_DOMAIN,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Domain~Citation_Type_Label)

# English vs. Non-English (Delta - Core and Periphery) -------------------------

Delta_Degree_Centrality_df_All_CorePeriphery <- plyr::ldply(list.files(pattern="OUTPUT_Python_Field_Delta_Degree_Centrality_KLD_Citation_Corpus_RAKE_and_GoogleAPI_All_*",
                                                                       full.names=T), 
                                                            read.csv)  %>%
  dplyr::mutate(Country = stringr::str_to_title(stringr::str_replace(as.character(Country),"_"," "))) %>%
  subset(Percentile_Cutoff==0)%>%
  subset(Citation_Window==5)%>%
  subset(Citation_Type=="Deflated")%>%
  dplyr::select(Country,Year,Discipline,OutDelta_Citation_ji)%>%
  dplyr::rename(InDelta_Citation_ij = OutDelta_Citation_ji) %>%
  dplyr::mutate(Language = "All")%>%
  rbind(.,Delta_Degree_Centrality_df%>%
          subset(Citation_Type=="Field_Deflated")%>%
          subset(Citation_Window==5)%>%
          subset(Journal_Censoring=="Since_1980")%>%
          subset(Percentile_Cutoff==0)%>%
          dplyr::select(Country,Year,Discipline,InDelta_Citation_ij)%>%
          dplyr::mutate(Language = "English_Only"))%>%
  as.data.frame()%>%
  merge(.,subset(MAG_Field_ID_df,select=c("FieldID","Domain")),by.x=c("Discipline"),by.y=c("FieldID"),how="left")%>%
  subset(Year>=minimum_year & Year<=(maximum_year-Citation_Window_Criterion_)) %>% 
  dplyr::mutate(Language = case_when(as.character(Language)=="English_Only"~"English Only",
                                     TRUE~"English and\nTranslated")) %>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))

# English vs. Non-English (Delta - Core and Periphery) | Data and Plot -------------------------

Delta_Degree_Centrality_df_CorePeriphery_Language <- Delta_Degree_Centrality_df_All_CorePeriphery %>%
  dplyr::group_by(Year,Discipline,Language) %>%
  dplyr::mutate(InDelta_Citation_ij = mean_center(InDelta_Citation_ij)) %>%
  dplyr::mutate(Core_Periphery = case_when(Country %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core"))) %>%
  dplyr::group_by(Year,Core_Periphery,Language) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE = sd(InDelta_Citation_ij,na.rm=T)/sqrt(length(InDelta_Citation_ij)))

SM_Language_plot_Delta_over_CorePeriphery <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_Language ,
                                                    aes(x=Year,
                                                        y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  ### Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,
                  fill=Core_Periphery),
              alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_Language,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery), 
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_wrap(~Language)

# English vs. Non-English (Delta - Core and Periphery) | Domain Data and Plot -------------------------

Delta_Degree_Centrality_df_CorePeriphery_DOMAIN_Language <- Delta_Degree_Centrality_df_All_CorePeriphery %>%
  dplyr::group_by(Year,Discipline,Language,Domain) %>%
  dplyr::mutate(InDelta_Citation_ij = mean_center(InDelta_Citation_ij)) %>%
  dplyr::group_by(Year,Core_Periphery,Language,Domain) %>%
  dplyr::summarise(InDelta_Citation_ij_Mean=mean(InDelta_Citation_ij,na.rm = T),
                   InDelta_Citation_ij_SE = sd(InDelta_Citation_ij,na.rm=T)/sqrt(length(InDelta_Citation_ij)))%>%
  as.data.frame()

SM_Language_plot_Delta_over_CorePeriphery_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_DOMAIN_Language ,
                                                           aes(x=Year,
                                                               y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  ### Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,
                  fill=Core_Periphery),
              alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_DOMAIN_Language,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery), 
  #                           colour="white",fontface = 'bold',segment.color = "black")+
  
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Domain~Language)

# Journal Censored vs. All (Delta - Core and Periphery) -------------------------
# Journal Censored vs. All (Delta - Core and Periphery) | Data and Plot -------------------------

Delta_Degree_Centrality_df_CorePeriphery_JournalCensored <- Delta_Degree_Centrality_df_CorePeriphery_SM %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  #subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Journal_Censoring_Label = case_when(Journal_Censoring=="No_Censoring"~"No Journal\nCensoring",
                                                    TRUE~"Journals Censored\nSince 1980"))%>%
  as.data.frame()

SM_JournalCensored_plot_Delta_over_CorePeriphery <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_JournalCensored,
                                                           aes(x=Year,
                                                               y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,fill=Core_Periphery),alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_JournalCensored,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(~Journal_Censoring_Label)

# Journal Censored vs. All (Delta - Core and Periphery) | Domain Data and Plot -------------------------

Delta_Degree_Centrality_df_CorePeriphery_JournalCensored_DOMAIN <- Delta_Degree_Centrality_df_CorePeriphery_DOMAIN_SM %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  #subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_)%>%
  dplyr::mutate(Journal_Censoring_Label = case_when(Journal_Censoring=="No_Censoring"~"No Journal\nCensoring",
                                                    TRUE~"Journals Censored\nSince 1980"))%>%
  as.data.frame()

SM_JournalCensored_plot_Delta_over_CorePeriphery_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_CorePeriphery_JournalCensored_DOMAIN,
                                                                  aes(x=Year,
                                                                      y=InDelta_Citation_ij_Mean))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed") +
  #geom_vline(xintercept = 2000,colour="black",linetype="dotted")+
  #geom_vline(xintercept = maximum_year-Citation_Window_Criterion_,colour="black",linetype="dotted")+
  theme_bw() + 
  ylab("Average Citational Distortion\n(Standard Deviations)")+
  xlab("Years")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="bottom")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  
  # Example Countries
  geom_point(aes(colour=Core_Periphery))+
  geom_line(aes(colour=Core_Periphery))+
  geom_ribbon(aes(ymax=InDelta_Citation_ij_Mean+InDelta_Citation_ij_SE,
                  ymin=InDelta_Citation_ij_Mean-InDelta_Citation_ij_SE,
                  colour=Core_Periphery,fill=Core_Periphery),alpha=0.2)+
  
  # ggrepel::geom_label_repel(data=subset(Delta_Degree_Centrality_df_CorePeriphery_JournalCensored_DOMAIN,
  #                                       Year==2010),
  #                           aes(fill=Core_Periphery,label=Core_Periphery),
  #                           fontface = 'bold', color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Domain~Journal_Censoring_Label,scales="free_y")

# Plot 4 | Delta Centrality 2000 v 2012  ---------------------------------
# Plot 4 | Delta Centrality 2000 v 2012 | Data and Plot ---------------------------------


Delta_Degree_Centrality_df_Trends_20002012 <- Delta_Degree_Centrality_df_Trends %>%
  dplyr::mutate(Percentile_Cutoff_Labels = case_when(Percentile_Cutoff==0 ~ "All Inclusive",
                                                     Percentile_Cutoff==0.25 ~ "Above 25th\nPercentile",
                                                     Percentile_Cutoff==0.75 ~ "Above 75th\nPercentile")) %>%
  #subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(Citation_Type==Inflation_Criterion_) %>%
  subset(Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(Year>=minimum_year & Year<=maximum_year-Citation_Window_Criterion_) %>%
  dplyr::mutate(Country = tools::toTitleCase(tolower(stringr::str_replace(Country,"_"," ")))) %>%
  ungroup() %>%
  dplyr::mutate(Nation_Name = case_when(Country == "United States" ~ "USA",
                                        Country == "United Kingdom" ~ "UK",
                                        Country == "South Korea" ~ "South Korea",
                                        Country == "South Africa" ~ "South Africa",
                                        Country == "Democratic Republic Of_the_congo" ~ "DR Congo",
                                        TRUE ~ stringr::str_to_title(stringr::str_replace(as.character(Country),"_"," ")))) %>%
  as.data.frame() %>%
  subset(select=c(Nation_Name,Year,InDelta_Citation_ij_Mean,Percentile_Cutoff)) %>%
  subset(Year%in%c(2000,maximum_year-Citation_Window_Criterion_)) %>%
  tidyr::spread(Year,InDelta_Citation_ij_Mean) %>%
  tidyr::drop_na()

# Plot 4 | Delta Centrality 2000 v 2012 | Main and SM Cutoff | Plot  ---------------------------------

plot_Delta_2000_2012_list <- list()
i <- 1
for(cutoff_ in c(0,0.25,0.75)){
  
  Delta_Trend_Data <- Delta_Degree_Centrality_df_Trends_20002012%>%
    subset(Percentile_Cutoff==cutoff_)%>%
    subset(select=-c(Percentile_Cutoff)) %>%
    dplyr::mutate(Combined=`2000`+`2012`)
  
  Delta_Trend_Data_Top <- Delta_Trend_Data %>% 
    top_n(10,Combined) %>%
    subset(`2000`>0 & `2012`>0) %>%
    mutate(QuadColor = "Upper")
  
  Delta_Trend_Data_Bottom <- Delta_Trend_Data %>% 
    top_n(-10,Combined) %>%
    subset(`2000`<0 & `2012`<0) %>%
    mutate(QuadColor = "Lower")
  
  Delta_Trend_Highlights <- rbind(Delta_Trend_Data_Top,Delta_Trend_Data_Bottom)
  Delta_Trend_Highlights$QuadColor <- as.factor(Delta_Trend_Highlights$QuadColor)
  
  
  cutoff_label_ <- dplyr::case_when(cutoff_==0 ~ "All Inclusive",
                                    cutoff_==0.25 ~ "Above 25th\nPercentile",
                                    cutoff_==0.75 ~ "Above 75th\nPercentile")
  
  plot_Delta_2000_2012 <- ggplot(data=Delta_Trend_Data,
                                 aes(x=`2000`,y=`2012`,label=Nation_Name))+
    geom_point(alpha=0.2)+
    geom_hline(yintercept = 0,colour="red",linetype="dashed")+
    geom_vline(xintercept = 0,colour="red",linetype="dashed")+
    geom_abline(a=1)+
    #coord_trans(y = "log1p",x = "log1p")+
    theme_bw() + 
    ylab("Average Citational Distortion (SDs)\n2012")+
    xlab("Average Citational Distortion (SDs)\n2000")+
    theme(text = element_text(size=16),
          panel.border = element_blank(), 
          panel.grid.major = element_blank(),
          panel.grid.minor = element_blank(), 
          axis.line = element_line(colour = "black")) +
    theme(legend.position="None")+
    theme(strip.background = element_blank(), strip.placement = "outside")+
    theme(legend.title=element_blank())+
    theme(legend.background = element_rect(colour = "black"))+
    ggtitle(paste(as.character(cutoff_label_)))+
    theme(plot.title = element_text(hjust = 0.5))+
    #gghighlight::gghighlight((`2000`>Delta_Trend_Data_Top$top2000 & `2012`>Delta_Trend_Data_Top$top2012) | 
    #                        (`2000`<Delta_Trend_Data_Bottom$bottom2000 & `2012`< Delta_Trend_Data_Bottom$bottom2012),
    #                     label_key = Nation_Name)+
    
    geom_point(data=Delta_Trend_Highlights)+
    ggrepel::geom_label_repel(data=Delta_Trend_Highlights,
                              fontface="bold",
                              color="white",
                              #fill = "white",
                              show.legend = F,
                              aes(fill=QuadColor))+
    
    ggpubr::stat_cor(group=1,method = "pearson",label.x.npc = 'left', label.y.npc = "top",label.sep="\n")+
    ggsci::scale_fill_npg()
  
  plot_Delta_2000_2012_list[[i]] <- plot_Delta_2000_2012
  i <- i + 1
}


SM_Cutoff_plot_Delta_2000_2012 <- ggpubr::ggarrange(plot_Delta_2000_2012_list[[2]],
                                                    plot_Delta_2000_2012_list[[3]],
                                                    #plot_Delta_2000_2012_list[[4]],
                                                    #plot_Delta_2000_2012_list[[5]],
                                                    ncol=2,
                                                    nrow = 1)

plot_Delta_2000_2012_Percentile <- plot_Delta_2000_2012_list[[1]] + 
  ggtitle("")

# Plot 4A | Delta Centrality 2000 v 2012 Region | Data  and Plot ---------------------------------


Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION <- Delta_Degree_Centrality_df_Domain %>%
  subset(Citation_Window==Citation_Window_Criterion_ & 
           Citation_Type==Inflation_Criterion_&
           Journal_Censoring==Journal_Censoring_Criterion_) %>%
  subset(select=-c(Citation_Window,Citation_Type,Journal_Censoring)) %>%
  merge(.,Regional_df,
        by=c("Country"))%>%
  #dplyr::mutate(Country = tools::toTitleCase(tolower(stringr::str_replace(Country,"_"," ")))) %>%
  ungroup() %>%
  dplyr::mutate(Nation_Name = case_when(Country == "United States" ~ "USA",
                                        Country == "United Kingdom" ~ "UK",
                                        Country == "South Korea" ~ "South Korea",
                                        Country == "South Africa" ~ "South Africa",
                                        Country == "Democratic Republic_of_the_congo" ~ "DR Congo",
                                        TRUE ~ stringr::str_to_title(stringr::str_replace(as.character(Country),"_"," ")))) %>%
  as.data.frame() %>%
  subset(select=c(Percentile_Cutoff,Domain,Region,Nation_Name,Year,InDelta_Citation_ij_Mean)) %>%
  subset(Year%in%c(2000,maximum_year-Citation_Window_Criterion_)) %>%
  tidyr::spread(Year,InDelta_Citation_ij_Mean) %>%
  tidyr::drop_na() %>%
  subset(!(Region%in%c("US and Canada","Oceania")))


# Numbers for European Correlation Reference

Values_for_Europe_Correlation_2000_2012_df <- Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION %>%
  subset(Region=="Europe")%>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  select(-c(Percentile_Cutoff,Region,Nation_Name)) %>%
  dplyr::group_by(Domain)%>%
  dplyr::summarise(Rho=cor.test(`2000`,`2012`,method = "pearson")$estimate,
                   t = cor.test(`2000`,`2012`,method = "pearson")$statistic,
                DF = cor.test(`2000`,`2012`,method = "pearson")$parameter,
                CI_L = cor.test(`2000`,`2012`,method = "pearson")$conf.int[1],
                CI_H = cor.test(`2000`,`2012`,method = "pearson")$conf.int[2],
                PValue=cor.test(`2000`,`2012`,method = "pearson")$p.value) %>%
  as.data.frame()

Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION_Percentile <- Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  select(-c(Percentile_Cutoff)) 
  #dplyr::mutate(Nation_Name_and_Coordinate = paste(Nation_Name," (",round(`2000`,2),",",round(`2012`,2),")",sep=""))

plot_Delta_2000_2012_DOMAIN_REGION <- ggplot(data=Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION_Percentile,
                                             aes(x=`2000`,y=`2012`))+
  geom_point(aes(colour=Region))+
  geom_hline(yintercept = 0,colour="red",linetype="dashed")+
  geom_vline(xintercept = 0,colour="red",linetype="dashed")+
  geom_abline(slope=1)+
  #coord_trans(y = "log1p",x = "log1p")+
  theme_bw() + 
  ylab("Average Citational Distortion (SDs)\n2012")+
  xlab("Average Citational Distortion (SDs)\n2000")+
  theme(text = element_text(size=16),
        panel.border = element_blank(), 
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(), 
        axis.line = element_line(colour = "black")) +
  theme(legend.position="none")+
  theme(strip.background = element_blank(), strip.placement = "outside")+
  theme(legend.title=element_blank())+
  theme(legend.background = element_rect(colour = "black"))+
  ggtitle("")+
  theme(plot.title = element_text(hjust = 0.5))+
  ggsci::scale_color_npg()+
  #gghighlight::gghighlight( 
  #                           (`2000`< 0.1 & `2012`> 0.1 & Region=="Asia"),
  #                         label_key = Nation_Name)+
  ggpubr::stat_cor(group=1,method = "pearson",label.x.npc = 'left', label.y.npc = "top",label.sep='\n')+
   # ggrepel::geom_label_repel(data=Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION_Percentile%>%
   #                    subset(Nation_Name%in%c("China","India","France","Germany","Mexico","Brazil","South Africa","Botswana")),
   #                  aes(fill=Region,label=Label), 
   #                  size =2,
   #                  direction="both",
   #                  fontface = 'bold', 
   #                  force = 5,
   #                  color = 'white',segment.color = "black")+
  ggsci::scale_fill_npg()+
  facet_grid(Domain~Region)

# Plot 4B | Delta Centrality 2000 v 2012 Core v Periphery Percent Aveages | Plot -------------------------------------

CorePeriphey_Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION_Percentile <- Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION_Percentile %>%
  tidyr::pivot_longer(!c("Domain","Region","Nation_Name"), names_to = "Year", values_to = "Distortion") %>%
  
  dplyr::mutate(Core_Periphery = case_when(Nation_Name %in% Censored_Core_Country_List~"Core",
                                           TRUE ~ "Periphery"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Core","Periphery")))%>%
  
  dplyr::group_by(Domain,Year,Core_Periphery)%>%
  dplyr::mutate(Over_Undercited = case_when(Distortion>0  ~ "Overcited",
                                            Distortion<0  ~ "Undercited",
                                            TRUE ~ as.character("Error")))%>%
  dplyr::group_by(Domain,Year,Core_Periphery,Over_Undercited)%>%
  dplyr::summarise(N = length(Distortion),
                   Distortion_Mean = mean(Distortion,na.rm=T),
                   Distortion_SE = sd(Distortion,na.rm=T)/sqrt(length(Distortion)))%>%
  dplyr::group_by(Domain,Core_Periphery)%>%
  dplyr::mutate(Percent = N/sum(N)) %>%
  dplyr::mutate(Labels = paste(paste(as.character(round(Percent*100,2)),"% ",paste("N=",as.character(N),sep=""),sep=""),
                               paste(as.character(round(Distortion_Mean,2))," (",as.character(round(Distortion_SE,2)),")",sep=""),
                               sep="\n"))%>%
  dplyr::mutate(Core_Periphery = factor(Core_Periphery,levels=c("Periphery","Core")))%>% #swapped so that red is periphery
  as.data.frame()

make_labels <- function(labels) {
  result <- stringr::str_replace(labels, "\\.","\\\n")
  return(result)
}

plot_Delta_2000_2012_Core_Periphery_DOMAIN_REGION <- ggplot(data=CorePeriphey_Delta_Degree_Centrality_Trends_df_20002012_DOMAIN_REGION_Percentile,
       aes(x=interaction(Year),y=Percent))+
  geom_line(aes(group=interaction(Core_Periphery),colour=Core_Periphery),linetype="dashed")+
  geom_point(aes(group=interaction(Core_Periphery),colour=Core_Periphery),show.legend = F)+
  ggrepel::geom_label_repel(aes(label=Labels,
                                fill=Core_Periphery),
                   colour="white",
                   #nudge_y = 0.15,
                   #nudge_x=0.25,
                   fontface="bold",
                   segment.color="black",
                   show.legend = F)+
  theme_bw() + 
     ylab("Percent of Countries")+
     xlab("")+
     theme(text = element_text(size=16),
                     panel.border = element_blank(), 
                     panel.grid.major = element_blank(),
                     panel.grid.minor = element_blank(), 
                     axis.line = element_line(colour = "black")) +
     theme(legend.position="bottom")+
     theme(strip.background = element_blank(), strip.placement = "outside")+
     theme(legend.title=element_blank())+
     theme(legend.background = element_rect(colour = "black"))+
     ggtitle("")+
     theme(plot.title = element_text(hjust = 0.5))+
  scale_x_discrete(labels = make_labels, name = "")+
  scale_y_continuous(labels = scales::percent)+
  guides(color = guide_legend(override.aes = list(size = 10)))+
  ggsci::scale_fill_npg()+
  ggsci::scale_color_npg()+
  facet_grid(Over_Undercited~Domain,scale="free_y")

# Plot 4C | Delta Centrality 2000 v 2012 World Map | Data  ---------------------------------

require(maps)
require(viridis)
world_map <- map_data("world")

Delta_Degree_Centrality_df_Map_Region_Domain <- Delta_Degree_Centrality_df_Domain %>%
  subset(Citation_Window==Citation_Window_Criterion_ & 
           Citation_Type==Inflation_Criterion_ & 
           Year>=minimum_year & 
           Year<=maximum_year-Citation_Window_Criterion_) %>%
  subset(Percentile_Cutoff==Percentile_Cutoff_Criterion_) %>%
  subset(select=-c(Percentile_Cutoff)) %>%
  merge(Regional_df,.,by=c("Country"),all.x=T) %>%
  dplyr::mutate(Country = tools::toTitleCase(tolower(stringr::str_replace(Country,"_"," ")))) %>%
  ungroup() %>%
  dplyr::mutate(Nation_Name = case_when(Country == "United States" ~ "USA",
                                        Country == "United Kingdom" ~ "UK",
                                        Country == "South Korea" ~ "South Korea",
                                        Country == "South Africa" ~ "South Africa",
                                        Country == "Democratic Republic_of_the_congo" ~ "DR Congo",
                                        TRUE ~ stringr::str_to_title(stringr::str_replace(as.character(Country),"_"," ")))) %>%
  as.data.frame() %>%
  subset(select=c(Domain,Nation_Name,Year,InDelta_Citation_ij_Mean,Region))


total_world_map = rbind(
  data.frame(world_map,Year=maximum_year-Citation_Window_Criterion_))

# Plot 4C | Delta Centrality 2000 v 2012 World Regional | Data  ---------------------------------

Delta_Degree_Centrality_df_Map_Coord_Region_Domain <- Delta_Degree_Centrality_df_Map_Region_Domain %>%
  subset(Year%in%c(maximum_year-Citation_Window_Criterion_)) %>%
  subset(Region %in% c("Europe","MENA and\nSub-Saharan Africa","Latin America\nand Caribbean","Asia")) %>%
  dplyr::group_by(Domain,Region,Nation_Name,Year) %>%
  dplyr::summarise(InDelta_Citation_ij = mean(InDelta_Citation_ij_Mean,na.rm = T)) %>%
  merge(total_world_map,.,
        by.x=c("region","Year"),
        by.y=c("Nation_Name","Year"),
        all.x=T) %>%
  subset(region!="Antarctica") %>%
  dplyr::arrange(group,order) 

# Plot 4C | Global | Plot ---------

Delta_Degree_Centrality_df_Map_Coord_Region_Domain_GLOBE <- Delta_Degree_Centrality_df_Map_Region_Domain %>%
  subset(Year%in%c(maximum_year-Citation_Window_Criterion_)) %>%
  dplyr::group_by(Domain,Region,Nation_Name,Year) %>%
  dplyr::summarise(InDelta_Citation_ij = mean(InDelta_Citation_ij_Mean,na.rm = T)) %>%
  merge(total_world_map,.,
        by.x=c("region","Year"),
        by.y=c("Nation_Name","Year"),
        all.x=T) %>%
  subset(region!="Antarctica") %>%
  dplyr::arrange(group,order) %>%
  dplyr::mutate(Region = case_when(region == "Bosnia and Herzegovina" ~ "Europe",
                                   region == "Czech Republic" ~ "Europe",
                                   region == "South Sudan" ~ "Africa",
                                   region == "Kyrgyzstan" ~ "Asia",
                                   region == "Central African Republic" ~ "Africa",
                                   region == "Ivory Coast" ~ "Africa",
                                   region == "Democratic Republic of the Congo" ~ "Africa",
                                   region == "Republic of Congo" ~ "Africa",
                                   TRUE ~ as.character(Region))) %>%
  dplyr::mutate(Region = dplyr::case_when(region=="USA"~"North\nAmerica",
                                          region=="Canada"~"North\nAmerica",
                                          region=="Australia"~"Oceania",
                                          region=="New Zealand"~"Oceania",
                                          TRUE~Region))%>%
  dplyr::mutate(InDelta_Citation_ij = dplyr::case_when(
    is.na(InDelta_Citation_ij)==T~as.numeric(NA),
    TRUE~as.numeric(InDelta_Citation_ij))) %>%
  dplyr::group_by(Year,Region) %>%
  dplyr::mutate(InDelta_Citation_ij_Region = mean(InDelta_Citation_ij,na.rm=T)) %>%
  subset(is.na(InDelta_Citation_ij_Region)==F)%>%
  subset(is.na(Region)==F)%>%
  subset(is.na(Domain)==F)

plot_Global_DeltaCentrality_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_Map_Coord_Region_Domain_GLOBE, 
       aes(map_id = region, 
           fill = InDelta_Citation_ij_Region))+
  geom_map(map = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_GLOBE,  color = "white",size=0)+
  expand_limits(x = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_GLOBE$long, 
                y = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_GLOBE$lat)+
  scale_fill_viridis_c(option = "C",direction=-1)+
  theme_void()+
  theme(legend.position="left")+  
  theme(strip.text.x = element_text(size = 10))+
  theme(strip.text.y = element_text(size = 10))+
  labs(fill='Average Citational\nDistortion (SDs)')+
  facet_wrap(~Domain,nrow=2)

# Plot 4D | Delta Centrality 2000 v 2012 World Regional | Europe | Plot  ---------------------------------

Delta_Degree_Centrality_df_Map_Coord_Region_Domain_EUROPE <- Delta_Degree_Centrality_df_Map_Coord_Region_Domain %>% 
  subset(Region%in%c("Europe")) %>%
  dplyr::mutate(InDelta_Citation_ij = dplyr::case_when(
    region=="Russia"~as.numeric(NA),
    TRUE~InDelta_Citation_ij)) %>%
  subset(!region%in%c("Uzbekistan","Tajikistan","Kazakhstan","Azerbaijan","Greenland"))

plot_Europe_DeltaCentrality_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_Map_Coord_Region_Domain_EUROPE, 
                                             aes(map_id = region, 
                                                 fill = InDelta_Citation_ij))+
  geom_map(map = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_EUROPE,  color = "white",size=0)+
  expand_limits(x = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_EUROPE$long, 
                y = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_EUROPE$lat)+
  scale_fill_viridis_c(option = "C",direction=-1)+
  theme_void()+
  theme(legend.position="left")+  
  theme(strip.text.x = element_text(size = 10))+
  theme(strip.text.y = element_text(size = 10))+
  labs(fill='Average Citational\nDistortion (SDs)')+
  facet_grid(~Domain)

# Plot 4E | Delta Centrality 2000 v 2012 World Regional | Asia | Plot  ---------------------------------


Delta_Degree_Centrality_df_Map_Coord_Region_Domain_ASIA <- Delta_Degree_Centrality_df_Map_Coord_Region_Domain %>% 
  subset(Region%in%c("Asia")) %>%
  as.data.frame()

plot_Asia_DeltaCentrality_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_Map_Coord_Region_Domain_ASIA, 
                                           aes(map_id = region, 
                                               fill = InDelta_Citation_ij))+
  geom_map(map = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_ASIA,  color = "white",size=0)+
  expand_limits(x = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_ASIA$long, 
                y = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_ASIA$lat)+
  scale_fill_viridis_c(option = "C",direction=-1)+
  theme_void()+
  theme(legend.position="left")+  
  theme(strip.text.x = element_text(size = 10))+
  theme(strip.text.y = element_text(size = 10))+
  labs(fill='Average Citational\nDistortion (SDs)')+
  facet_grid(~Domain)

# Plot 4F | Delta Centrality 2000 v 2012 World Regional | Africa and MENA | Plot  ---------------------------------


Delta_Degree_Centrality_df_Map_Coord_Region_Domain_AFRICAMENA <- Delta_Degree_Centrality_df_Map_Coord_Region_Domain %>%
  subset(Region %in% c("MENA and\nSub-Saharan Africa")) %>%
  as.data.frame()


plot_AfricaMENA_DeltaCentrality_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_Map_Coord_Region_Domain_AFRICAMENA, 
                                                 aes(map_id = region, 
                                                     fill = InDelta_Citation_ij))+
  geom_map(map = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_AFRICAMENA,  
           color = "white",size=0)+
  expand_limits(x = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_AFRICAMENA$long, 
                y = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_AFRICAMENA$lat)+
  scale_fill_viridis_c(option = "C",direction=-1)+
  theme_void()+
  theme(legend.position="left")+  
  theme(strip.text.x = element_text(size = 10))+
  theme(strip.text.y = element_text(size = 10))+
  labs(fill='Average Citational\nDistortion (SDs)')+
  facet_grid(~Domain)

# Plot 4G | Delta Centrality 2000 v 2012 World Regional | Lat Am | Plot  ---------------------------------

Delta_Degree_Centrality_df_Map_Coord_Region_Domain_LATAM <- Delta_Degree_Centrality_df_Map_Coord_Region_Domain %>%
  subset(Region %in% c("Latin America\nand Caribbean")) %>%
  as.data.frame()

plot_LatAM_DeltaCentrality_DOMAIN <- ggplot(data=Delta_Degree_Centrality_df_Map_Coord_Region_Domain_LATAM, 
                                            aes(map_id = region, 
                                                fill = InDelta_Citation_ij))+
  geom_map(map = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_LATAM,  color = "white",size=0)+
  expand_limits(x = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_LATAM$long, 
                y = Delta_Degree_Centrality_df_Map_Coord_Region_Domain_LATAM$lat)+
  scale_fill_viridis_c(option = "C",direction=-1)+
  theme_void()+
  theme(legend.position="left")+  
  theme(strip.text.x = element_text(size = 10))+
  theme(strip.text.y = element_text(size = 10))+
  labs(fill='Average Citational\nDistortion (SDs)')+
  facet_grid(~Domain)


# Combine and Output Plots  ---------------------------------


plot_3 <- ggpubr::ggarrange(plot_Delta_over_CountrySpecific,
                              plot_Delta_CorePeriphery,
                              labels=c("A","B"),
                              nrow=1,
                              ncol=2) 

plot_4 <- ggpubr::ggarrange(plot_Delta_over_CountrySpecific_Domain,
                             plot_Delta_CorePeriphery_Domain,
                             labels=c("A","B"),
                             nrow=2,
                             ncol=1) 

# NOTE: NO SM
plot_5 <- ggpubr::ggarrange(plot_Delta_2000_2012_DOMAIN_REGION,
                             nrow=1,
                             ncol=1) 

plot_6 <- ggpubr::ggarrange(plot_Delta_2000_2012_Core_Periphery_DOMAIN_REGION,
                             nrow=1,
                             ncol=1) 

plot_7 <- ggpubr::ggarrange(plot_Global_DeltaCentrality_DOMAIN,
                             nrow=1,
                             ncol=1) 

plot_8 <- ggpubr::ggarrange(plot_Europe_DeltaCentrality_DOMAIN,
                              plot_Asia_DeltaCentrality_DOMAIN,
                              plot_LatAM_DeltaCentrality_DOMAIN,
                              plot_AfricaMENA_DeltaCentrality_DOMAIN,
                              labels=c("A","B","C","D"),
                              nrow=4,
                              ncol=1)

plot_8ab <- ggpubr::ggarrange(plot_Europe_DeltaCentrality_DOMAIN,
                              plot_Asia_DeltaCentrality_DOMAIN,
                              labels=c("A","B"),
                              nrow=2,
                              ncol=1) 

plot_8cd <- ggpubr::ggarrange(plot_LatAM_DeltaCentrality_DOMAIN,
                              plot_AfricaMENA_DeltaCentrality_DOMAIN,
                              labels=c("C","D"),
                              nrow=2,
                              ncol=1) 



pdf(paste("NEW_FIGURE_R_Plot3.pdf",sep=""),width=(8*1.75), height=8*1)
plot(plot_3)
dev.off()

pdf(paste("NEW_FIGURE_R_Plot4.pdf",sep=""),width=(10*1.25), height=10*1)
plot(plot_4)
dev.off()

pdf(paste("NEW_FIGURE_R_Plot5.pdf",sep=""),width=(8*1.25), height=8*1.25)
plot(plot_5)
dev.off()

pdf(paste("NEW_FIGURE_R_Plot6.pdf",sep=""),width=(8*1.65), height=8*1)
plot(plot_6)
dev.off()

pdf(paste("NEW_FIGURE_R_Plot7.pdf",sep=""),width=(8*1.5), height=8*1.)
plot(plot_7)
dev.off()


pdf(paste("NEW_FIGURE_R_Plot8.pdf",sep=""),width=(8*1.5), height=8*2)
plot(plot_8)
dev.off()


pdf(paste("SM_Figure_Coherence_Cutoffs_Delta_2000_2012.pdf",sep=""),width=(8*1.5), height=8*1)
plot(SM_Cutoff_plot_Delta_2000_2012)
dev.off()

pdf(paste("SM_Figure_Coherence_Cutoffs_Delta_over_CountrySpecific.pdf",sep=""),width=(8*1.75), height=8*1)
plot(SM_Cutoff_plot_Delta_over_CountrySpecific)
dev.off()

pdf(paste("SM_Figure_Coherence_Cutoffs_Delta_over_CountrySpecific_DOMAIN.pdf",sep=""),width=(8*1.5), height=8*1.0)
plot(SM_Cutoff_plot_Delta_over_CountrySpecific_Domain)
dev.off()

pdf(paste("SM_Figure_Deflated_Delta_over_CountrySpecific.pdf",sep=""),width=(8*1.75), height=8*1.0)
plot(SM_Deflated_plot_Delta_over_CountrySpecific)
dev.off()

pdf(paste("SM_Figure_Deflated_Delta_over_CountrySpecific_DOMAIN.pdf",sep=""),width=(8*1.25), height=8*1.25)
plot(SM_Deflated_plot_Delta_over_CountrySpecific_DOMAIN)
dev.off()

pdf(paste("SM_Figure_JournalCensored_over_CountrySpecific.pdf",sep=""),width=(8*1.75), height=8*1)
plot(SM_JournalCensored_plot_Delta_over_CountrySpecific)
dev.off()

pdf(paste("SM_Figure_JournalCensored_over_CountrySpecific_Domain.pdf",sep=""),width=(8*1.25), height=8*1.25)
plot(SM_JournalCensored_plot_Delta_over_CountrySpecific_DOMAIN)
dev.off()

pdf(paste("SM_Figure_Language_over_CountrySpecific.pdf",sep=""),width=(8*1.75), height=8*1)
plot(SM_Language_plot_Delta_over_CountrySpecific)
dev.off()

pdf(paste("SM_Figure_Language_over_CountrySpecific_Domain.pdf",sep=""),width=(8*1.25), height=8*1.25)
plot(SM_Language_plot_Delta_over_CountrySpecific_DOMAIN)
dev.off()

pdf(paste("SM_Figure_Coherence_Cutoffs_Delta_over_CorePeriphery.pdf",sep=""),width=(8*1.75), height=8*1)
plot(SM_Cutoff_plot_Delta_over_CorePeriphery)
dev.off()

pdf(paste("SM_Figure_Coherence_Cutoffs_Delta_over_CorePeriphery_DOMAIN.pdf",sep=""),width=(8*1.5), height=8*1.0)
plot(SM_Cutoff_plot_Delta_over_CorePeriphery_DOMAIN)
dev.off()

pdf(paste("SM_Figure_Deflated_Delta_over_CorePeriphery.pdf",sep=""),width=(8*1.75), height=8*1.0)
plot(SM_Deflated_plot_Delta_over_CorePeriphery)
dev.off()

pdf(paste("SM_Figure_Deflated_Delta_over_CorePeriphery_DOMAIN.pdf",sep=""),width=(8*1.25), height=8*1.25)
plot(SM_Deflated_plot_Delta_over_CorePeriphery_DOMAIN)
dev.off()

pdf(paste("SM_Figure_Language_over_CorePeriphery.pdf",sep=""),width=(8*1.75), height=8*1)
plot(SM_Language_plot_Delta_over_CorePeriphery)
dev.off()

pdf(paste("SM_Figure_Language_over_CorePeriphery_Domain.pdf",sep=""),width=(8*1.25), height=8*1.25)
plot(SM_Language_plot_Delta_over_CorePeriphery_DOMAIN)
dev.off()

pdf(paste("SM_Figure_JournalCensored_over_CorePeriphery.pdf",sep=""),width=(8*1.75), height=8*1)
plot(SM_JournalCensored_plot_Delta_over_CorePeriphery)
dev.off()

pdf(paste("SM_Figure_JournalCensored_over_CorePeriphery_Domain.pdf",sep=""),width=(8*1.25), height=8*1.25)
plot(SM_JournalCensored_plot_Delta_over_CorePeriphery_DOMAIN)
dev.off()

pdf(paste("SM_Figure_FirstApperance_Impact_on_Delta_CorePeriphery.pdf",sep=""),width=(8*1.25), height=8*1)
plot(SM_plot_CorePeriphery_FirstApperance_Impact_on_Delta)
dev.off()

pdf(paste("SM_Figure_FirstApperance_Impact_on_Delta_CorePeriphery_Domain.pdf",sep=""),width=(10*1.25), height=10*1)
plot(SM_plot_CorePeriphery_FirstApperance_Impact_on_Delta_Domain)
dev.off()
